﻿#region Copyright information
// 
// Copyright © 2005-2013 Yongkee Cho. All rights reserved.
// 
// This code is a part of the Biological Object Library and governed under the terms of the
// GNU Lesser General  Public License (LGPL) version 2.1 which accompanies this distribution.
// For more information on the LGPL, please visit http://bol.codeplex.com/license.
// 
// - Filename: ExponentialDistribution.cs
// - Author: Yongkee Cho
// - Email: yongkeecho@gmail.com
// - Date Created: 2012-09-06 11:39 AM
// - Last Modified: 2013-01-25 3:59 PM
// 
#endregion
using System;
using System.Collections.Generic;
using System.Linq;
using BOL.Algorithms.Optimization;

namespace BOL.Maths.Distributions
{
    /// <summary>
    /// Represents an exponential distributin with scale inversely proportional to parameter lambda.
    /// <see>Mathematica: ref/ExponentialDistribution</see>
    /// <see>http://en.wikipedia.org/wiki/Exponential_distribution</see>
    /// </summary>
    public sealed class ExponentialDistribution : IUnivariateDistribution<double>, IEquatable<ExponentialDistribution>
    {
        #region Private variables

        private double _lambda;

        #endregion

        #region Public properties

        /// <summary>Gets or sets the lambda of an exponential distribution.</summary>
        public double Lambda
        {
            get { return _lambda; }
            set
            {
                if (value <= 0)
                    throw new ArgumentOutOfRangeException("value", "lambda is expected to be a positive real number.");
                
                _lambda = value;
            }
        }

        /// <summary>Gets the range of an exponential distribution.</summary>
        public IRange<double> Domain { get { return new Range<double>(0, double.PositiveInfinity); } }

        /// <summary>Gets the mean of an exponential distribution.</summary>
        public double Mean { get { return 1 / _lambda; } }

        /// <summary>Gets the median of an exponential distribution.</summary>
        public double Median { get { return Math.Log(2) / _lambda; } }

        /// <summary>Gets the mode of an exponential distribution.</summary>
        public double Mode { get { return 0; } }

        /// <summary>Gets the variance of an exponential distribution.</summary>
        public double Variance { get { return 1 / Math.Pow(_lambda, 2); } }

        /// <summary>Gets the skewness of an exponential distribution.</summary>
        public double Skewness { get { return 2; } }

        /// <summary>Gets the kurtosis of an exponential distribution.</summary>
        public double Kurtosis { get { return 6; } }

        /// <summary>Gets the entropy of an exponential distribution.</summary>
        public double Entropy { get { return 1 - Math.Log(_lambda); } }

        #endregion

        #region Constructor

        /// <summary>
        /// Instantiates an exponential distributin.
        /// </summary>
        public ExponentialDistribution(double lambda)
        {
            if (lambda <= 0)
                throw new ArgumentOutOfRangeException("lambda", "lambda is expected to be a positive real number.");
            
            _lambda = lambda;
        }

        #endregion

        #region ICloneable implementation

        public ExponentialDistribution Clone()
        {
            return new ExponentialDistribution(_lambda);
        }

        object ICloneable.Clone()
        {
            return Clone();
        }
        
        #endregion

        #region Public methods

        /// <summary>
        /// Returns the exponential density given value.
        /// </summary>
        /// <param name="value"></param>
        /// <returns>exponential density</returns>
        public double Pdf(double value)
        {
            if (value >= 0)
                return _lambda * Math.Exp(-_lambda * value);
            return 0;
        }

        /// <summary>
        /// Returns the cummulative exponential density given value.
        /// </summary>
        /// <param name="value"></param>
        /// <returns>cummulative exponential density</returns>
        public double Cdf(double value)
        {
            if (value >= 0)
                return 1 - Math.Exp(-_lambda * value);
            return 0;
        }

        /// <summary>
        /// Returns the value of the exponential distribution for given probability.
        /// </summary>
        /// <see cref="http://en.wikipedia.org/wiki/Quantile_function#Simple_example"/>
        /// <param name="p">probability</param>
        /// <returns></returns>
        public double Quantile(double p)
        {
            return -Math.Log(1 - p) / _lambda;
        }

        /// <summary>
        /// Returns a random value from the exponential distribution.
        /// <see>http://people.csail.mit.edu/milch/blog/index.html</see>
        /// <see>Reference: A Guide to Simulation, 2nd Ed. Bratley, Paul, Bennett L. Fox and Linus E. Schrage.</see>
        /// </summary>
        /// <returns>random value</returns>
        public double Sample(Random r)
        {
            return -Math.Log(r.NextDouble()) / _lambda;
        }

        /// <summary>
        /// Estimates parameters of the exponential distribution using maximum likelihood.
        /// <see>http://en.wikipedia.org/wiki/Exponential_distribution#Maximum_likelihood</see>
        /// </summary>
        public void MaximumLikelihoodEstimate(IEnumerable<double> source)
        {
            //ExponentialPosterior ep = new ExponentialPosterior(data, alpha, beta);
            //var optimizer = new IOptimizerByValue(ep, 5000, 1e-5, true, 0.1);
            //optimizer.Optimize();
            //Lambda = ep.Parameters[0];

            throw new NotImplementedException();
        }

        /// <summary>
        /// Estimates parameters of the exponential distribution using Bayesian methods.
        /// <see>http://en.wikipedia.org/wiki/Exponential_distribution#Bayesian_inference</see>
        /// </summary>
        public void BayesianEstimate(IEnumerable<double> source, params dynamic[] priorParameters)
        {
            throw new NotImplementedException();
        }

        #region Exponential posterior function

        private class ExponentialPosterior : ObjectiveFunctionBase
        {
            #region Private variables

            private readonly GammaDistribution _gamma; // a Gamma distribution

            #endregion

            #region Constructor

            public ExponentialPosterior(IEnumerable<double> data, double alpha, double beta)
                : base(1)
            {
                if (data == null)
                    throw new ArgumentNullException("data");

                var list = data.ToList();
                var length = list.Count;

                _gamma = new GammaDistribution(alpha + length, beta + length * list.Average());
            }

            #endregion

            #region Public methods

            public override double Value(double[] parameters)
            {
                return _gamma.Pdf(parameters[0]);
            }
            
            #endregion
        }

        #endregion

        #endregion

        #region IEquatable<ExponentialDistribution> implementation

        public bool Equals(ExponentialDistribution other)
        {
            return _lambda.Equals(other._lambda);
        }

        #endregion

        #region Object overriden

        public override int GetHashCode()
        {
            return Lambda.GetHashCode();
        }

        public override bool Equals(object other)
        {
            if (other == null)
                throw new ArgumentNullException();

            if (!(other is ExponentialDistribution))
                throw new InvalidCastException("The 'other' argument is not a ExponentialDistribution object.");
            
            return Equals(other as ExponentialDistribution);
        }

        public override string ToString()
        {
            return String.Format("Exponential(lambda = {0:G4})", _lambda);
        }

        #endregion
    }
}
